¿Qué es un modelo?

¿Qué es un modelo?

¿Qué es un modelo?

¿Qué es un modelo?

• Es habitual expresar un modelo de regresión lineal como

• Esto es idéntico a

• Por tanto, el modelo estadístico cuenta con:
  1. Una variable dependiente y una distribución de probabilidad,
  2. Los parámetros de la distribución de probabilidad seleccionada.
  3. Los parámetros estimados que, a veces, dependen de otras variables observables o latentes

Teoría clásica de los tests

La ecuación fundamental de la Psicometría:

• La puntuación verdadera es el valor esperado (i.e., la media) de la puntuación observada.

• Aunque \(X\), \(V\) y \(E\) son variables, la TCT no especifica ninguna distribución para ninguna de ellas.

Teoría clásica de los tests

La ecuación fundamental de la Psicometría:

\[ X = V + E \]

Teoría clásica de los tests

La ecuación fundamental de la Psicometría:

\[ X = V + E \]

• Se asume que \(E\sim\color{#7E57C2}{\mathcal{N}}(\color{#9a2515}\mu = 0,\color{#9a2515}\sigma = \color{#0197FD}{\sigma_E})\)

• Sin embargo, no se especifica la distribución de \(\color{#0197FD}{V}\).

Teoría moderna de los tests

Evolución natural de la TCT:

1. Modelos factoriales (Spearman 1904; Thurstone 1947; Jöreskog 1969)
2. Teoría de respuesta al Item (Lord, Novick, and Birnbaum 1968)

Teoría moderna de los tests

Evolución natural de la TCT:

1. Modelos factoriales (Spearman 1904; Thurstone 1947; Jöreskog 1969)
2. Teoría de respuesta al Item (Lord, Novick, and Birnbaum 1968)

Del valor esperado (\(\mathbb{E}[X]\)) a la causa latente (\(\eta\))

El modelo factorial confirmatorio

La respuesta de la persona \(i\) a un ítem se distribuye según

\[ \color{seagreen}{Y_i} \sim \color{#7E57C2}{\mathcal{N}}\left(\color{#9a2515}{\mu}=\color{#0197FD}{\beta_0} + \color{#0197FD}{\beta_1}\cdot \color{darkorange}{X_i} ,\,\color{#9a2515}{\sigma} = \color{#0197FD}{\sigma_\varepsilon} \right) \]

Teoría moderna de los tests

\[ \color{seagreen}{Y_i} \sim \color{#7E57C2}{\mathcal{N}}\left(\color{#9a2515}{\mu}=\color{#0197FD}{\beta_0} + \color{#0197FD}{\beta_1}\cdot \color{darkorange}{X_i} ,\,\color{#9a2515}{\sigma} = \color{#0197FD}{\sigma_\varepsilon} \right) \]

• Sin embargo, en psicometría:
  • se asume que la variable \(X\) no es observable; es una variable latente. Se suele representar con la letra griega \(\eta\).
  • suelen utilizarse las letras griegas \(\nu\) y \(\lambda\) en lugar de \(\beta_0\) y \(\beta_1\).
• Por eso, el mismo modelo se representa como

Teoría moderna de los tests

\[ \color{seagreen}{Y_i} \sim \color{#7E57C2}{\mathcal{N}}\left(\color{#9a2515}{\mu}=\color{#0197FD}{\beta_0} + \color{#0197FD}{\beta_1}\cdot \color{darkorange}{X_i} ,\,\color{#9a2515}{\sigma} = \color{#0197FD}{\sigma_\varepsilon} \right) \]

• Sin embargo, en psicometría:
  • se asume que la variable \(X\) no es observable; es una variable latente. Se suele representar con la letra griega \(\eta\).
  • suelen utilizarse las letras griegas \(\nu\) y \(\lambda\) en lugar de \(\beta_0\) y \(\beta_1\).
• Por eso, el mismo modelo se representa como

\[ \color{seagreen}{Y_i} \sim \color{#7E57C2}{\mathcal{N}}\left(\color{#9a2515}{\mu}=\color{#0197FD}{\nu} + \color{#0197FD}{\lambda} \cdot \color{darkorange}{\eta_i} ,\,\color{#9a2515}{\sigma} = \color{#0197FD}{\theta} \right) \]

Modelo factorial confirmatorio

\[ \color{seagreen}{Y_i} \sim \color{#7E57C2}{\mathcal{N}}\left(\color{#9a2515}{\mu}=\color{#0197FD}{\nu} + \color{#0197FD}{\lambda} \cdot \color{darkorange}{\eta_i} ,\,\color{#9a2515}{\sigma} = \color{#0197FD}{\theta} \right) \]

• \(\color{seagreen}{Y_i}\) es una medida imperfecta de una variable no observable, \(\color{darkorange}{\eta_i}\).
• \(\varepsilon_i\) recibe el nombre de unicidad (i.e., varianza no explicada).
  • Una parte de la unicidad es el error de medida
• Es habitual contar con varias medidas de la variable no observable.

Modelo factorial confirmatorio

\[ \color{seagreen}{Y_i} \sim \color{#7E57C2}{\mathcal{N}}\left(\color{#9a2515}{\mu}=\color{#0197FD}{\nu} + \color{#0197FD}{\lambda} \cdot \color{darkorange}{\eta_i} ,\,\color{#9a2515}{\sigma} = \color{#0197FD}{\theta} \right) \]

• \(\color{seagreen}{Y_i}\) es una medida imperfecta de una variable no observable, \(\color{darkorange}{\eta_i}\).
• \(\varepsilon_i\) recibe el nombre de unicidad (i.e., varianza no explicada).
  • Una parte de la unicidad es el error de medida
• Es habitual contar con varias medidas de la variable no observable.

\[ \begin{aligned} & Y_{i1} = \nu_1 + \lambda_1 \cdot \eta_i + \varepsilon_{i1},\quad\text{donde}\quad\varepsilon_{i1}\sim\mathcal{N}\left(\mu =0,\, \sigma = \theta_1\right)\\ & Y_{i2} = \nu_2 + \lambda_2 \cdot \eta_i + \varepsilon_{i2},\quad\text{donde}\quad\varepsilon_{i2}\sim\mathcal{N}\left(\mu =0,\, \sigma = \theta_2\right) \end{aligned} \]

Modelo factorial confirmatorio

\[ \color{seagreen}{Y_i} \sim \color{#7E57C2}{\mathcal{N}}\left(\color{#9a2515}{\mu}=\color{#0197FD}{\nu} + \color{#0197FD}{\lambda} \cdot \color{darkorange}{\eta_i} ,\,\color{#9a2515}{\sigma} = \color{#0197FD}{\theta} \right) \]

• \(\color{seagreen}{Y_i}\) es una medida imperfecta de una variable no observable, \(\color{darkorange}{\eta_i}\).
• \(\varepsilon_i\) recibe el nombre de unicidad (i.e., varianza no explicada).
  • Una parte de la unicidad es el error de medida
• Es habitual contar con varias medidas de la variable no observable.

\[ \begin{aligned} & Y_{i1} = \nu_1 + \lambda_1 \cdot \eta_i + \varepsilon_{i1},\quad\text{donde}\quad\varepsilon_{i1}\sim\mathcal{N}\left(\mu =0,\, \sigma = \theta_1\right)\\ & Y_{i2} = \nu_2 + \lambda_2 \cdot \eta_i + \varepsilon_{i2},\quad\text{donde}\quad\varepsilon_{i2}\sim\mathcal{N}\left(\mu =0,\, \sigma = \theta_2\right)\\ & Y_{i3} = \nu_3 + \lambda_3 \cdot \eta_i + \varepsilon_{i3},\quad\text{donde}\quad\varepsilon_{i3}\sim\mathcal{N}\left(\mu =0,\, \sigma = \theta_3\right)\\ \end{aligned} \]

Modelo factorial confirmatorio

\[ \color{seagreen}{Y_i} \sim \color{#7E57C2}{\mathcal{N}}\left(\color{#9a2515}{\mu}=\color{#0197FD}{\nu} + \color{#0197FD}{\lambda} \cdot \color{darkorange}{\eta_i} ,\,\color{#9a2515}{\sigma} = \color{#0197FD}{\theta} \right) \]

• \(\color{seagreen}{Y_i}\) es una medida imperfecta de una variable no observable, \(\color{darkorange}{\eta_i}\).
• \(\varepsilon_i\) recibe el nombre de unicidad (i.e., varianza no explicada).
  • Una parte de la unicidad es el error de medida
• Es habitual contar con varias medidas de la variable no observable.

\[ \begin{aligned} & Y_{i1} = \nu_1 + \lambda_1 \cdot \color{darkorange}{\eta_i} + \varepsilon_{i1},\quad\text{donde}\quad\varepsilon_{i1}\sim\mathcal{N}\left(\mu =0,\, \sigma = \theta_1\right)\\ & Y_{i2} = \nu_2 + \lambda_2 \cdot \color{darkorange}{\eta_i} + \varepsilon_{i2},\quad\text{donde}\quad\varepsilon_{i2}\sim\mathcal{N}\left(\mu =0,\, \sigma = \theta_2\right)\\ & Y_{i3} = \nu_3 + \lambda_3 \cdot \color{darkorange}{\eta_i} + \varepsilon_{i3},\quad\text{donde}\quad\varepsilon_{i3}\sim\mathcal{N}\left(\mu =0,\, \sigma = \theta_3\right)\\ \end{aligned} \]

Modelo factorial confirmatorio

\[ \color{seagreen}{Y_i} \sim \color{#7E57C2}{\mathcal{N}}\left(\color{#9a2515}{\mu}=\color{#0197FD}{\nu} + \color{#0197FD}{\lambda} \cdot \color{darkorange}{\eta_i} ,\,\color{#9a2515}{\sigma} = \color{#0197FD}{\theta} \right) \]

• \(\color{seagreen}{Y_i}\) es una medida imperfecta de una variable no observable, \(\color{darkorange}{\eta_i}\).
• \(\varepsilon_i\) recibe el nombre de unicidad (i.e., varianza no explicada).
  • Una parte de la unicidad es el error de medida
• Es habitual contar con varias medidas de la variable no observable.

Modelo factorial confirmatorio

\[ \color{seagreen}{Y_{i\color{black}{j}}} \sim \color{#7E57C2}{\mathcal{N}}\left(\color{#9a2515}{\mu}=\color{#0197FD}{\nu_\color{black}j} + \color{#0197FD}{\lambda_\color{black}j} \cdot \color{darkorange}{\eta_i} ,\,\color{#9a2515}{\sigma} = \color{#0197FD}{\theta_\color{black}j} \right) \]

• \(\color{seagreen}{Y_i}\) es una medida imperfecta de una variable no observable, \(\color{darkorange}{\eta_i}\).
• \(\varepsilon_i\) recibe el nombre de unicidad (i.e., varianza no explicada).
  • Una parte de la unicidad es el error de medida
• Es habitual contar con varias medidas de la variable no observable.

\[ Y_{i\color{red}j} = \nu_\color{red}j + \lambda_\color{red}j \cdot \eta_i + \varepsilon_{i\color{red}j},\quad\text{donde}\quad\varepsilon_{i\color{red}j}\sim\mathcal{N}\left(\mu =0,\, \sigma = \theta_\color{red}j\right) \]

• ¿Cómo ajustamos un modelo de regresión con una variable invisible?

• Respuesta: sorprendentemente, no la necesitamos.

Modelo factorial confirmatorio

• Los psicómetras prefieren la versión matricial del modelo.

• Hay una versión alternativa del modelo estadístico utilizando estas matrices.
• Esa versión no necesita las variables latentes

Modelo factorial confirmatorio

Finalmente, el modelo estadístico del factor común es

Modelo factorial confirmatorio

Finalmente, el modelo estadístico del factor común es

\[ \mathbf{y}_i \sim \mathcal{MVN}\left(\boldsymbol\mu = \color{#0197FD}{\boldsymbol\nu} + \color{#0197FD}{\mathbf{\Lambda}} \cdot \boldsymbol{\alpha},\,\mathbf{\Sigma} = \mathbf{\color{#0197FD}{\Lambda}}\cdot\mathbf{\Psi}\cdot\color{#0197FD}{\mathbf{\Lambda^\prime}} + \mathbf{\Theta}\right) \]

Modelo factorial confirmatorio

Finalmente, el modelo estadístico del factor común es

\[ \mathbf{y}_i \sim \mathcal{MVN}\left(\boldsymbol\mu = \boldsymbol\nu + \mathbf{\Lambda} \cdot \color{royalblue}{\boldsymbol{\alpha}},\,\mathbf{\Sigma} = \mathbf{\Lambda\cdot\color{royalblue}{\mathbf\Psi}\cdot\mathbf{\Lambda}^\prime} + \mathbf{\color{royalblue}{\mathbf{\Theta}}}\right) \]

• \(\color{royalblue}{\boldsymbol{\alpha}}\) es el vector con las medias de \(\boldsymbol{\eta}_i\)
• \(\color{royalblue}{\mathbf\Psi}\) es la matriz de covarianzas de \(\boldsymbol{\eta}_i\)
• \(\color{royalblue}{\mathbf\Theta}\) es la matriz de covarianzas de \(\boldsymbol{\varepsilon}_i\)

Modelo factorial confirmatorio

Finalmente, el modelo estadístico del factor común es

\[ \color{seagreen}{\mathbf{y}_i} \sim \color{#7E57C2}{\mathcal{MVN}}\left(\color{#9a2515}{\boldsymbol\mu} = \color{#0197FD}{\boldsymbol\nu} + \color{#0197FD}{\mathbf{\Lambda}} \cdot \color{#0197FD}{\boldsymbol{\alpha}},\,\color{#9a2515}{\mathbf{\Sigma}} = \color{#0197FD}{\mathbf{\Lambda}}\cdot\color{#0197FD}{\mathbf\Psi}\cdot\color{#0197FD}{\mathbf{\Lambda}^\prime} + \color{#0197FD}{\mathbf{\Theta}}\right) \]

• \(\boldsymbol{\alpha}\) es el vector con las medias de \(\boldsymbol{\eta}_i\)
• \(\mathbf\Psi\) es la matriz de covarianzas de \(\boldsymbol{\eta}_i\)
• \(\mathbf\Theta\) es la matriz de covarianzas de \(\boldsymbol{\varepsilon}_i\)

• \(\color{#0197FD}{\boldsymbol\nu} + \color{#0197FD}{\mathbf{\Lambda}} \cdot \color{#0197FD}{\boldsymbol{\alpha}}\) se denomina model-implied mean vector
• \(\color{#0197FD}{\mathbf{\Lambda}}\cdot\color{#0197FD}{\mathbf\Psi}\cdot\color{#0197FD}{\mathbf{\Lambda}^\prime} + \color{#0197FD}{\mathbf{\Theta}}\) se denomina model-implied covariance matrix

Modelo factorial confirmatorio

Finalmente, el modelo estadístico del factor común es

\[ \color{seagreen}{\mathbf{y}_i} \sim \color{#7E57C2}{\mathcal{MVN}}\left(\color{#9a2515}{\boldsymbol\mu} = \color{#0197FD}{\boldsymbol\nu} ,\,\color{#9a2515}{\mathbf{\Sigma}} = \color{#0197FD}{\mathbf{\Lambda}}\cdot\color{#0197FD}{\mathbf\Psi}\cdot\color{#0197FD}{\mathbf{\Lambda}^\prime} + \color{#0197FD}{\mathbf{\Theta}}\right) \]

• \(\boldsymbol{\alpha}\) es el vector con las medias de \(\boldsymbol{\eta}_i\)
• \(\mathbf\Psi\) es la matriz de covarianzas de \(\boldsymbol{\eta}_i\)
• \(\mathbf\Theta\) es la matriz de covarianzas de \(\boldsymbol{\varepsilon}_i\)

• \(\color{#0197FD}{\boldsymbol\nu} + \color{#0197FD}{\mathbf{\Lambda}} \cdot \color{#0197FD}{\boldsymbol{\alpha}}\) se denomina model-implied mean vector
• \(\color{#0197FD}{\mathbf{\Lambda}}\cdot\color{#0197FD}{\mathbf\Psi}\cdot\color{#0197FD}{\mathbf{\Lambda}^\prime} + \color{#0197FD}{\mathbf{\Theta}}\) se denomina model-implied covariance matrix

La edad según el frecuentismo

• La media de edad es una cantidad fija, pero desconocida (\(\mu_{\text{edad}}\)).
• El objetivo es estimar \(\mu_{\text{edad}}\) usando una muestra aleatoria.
  1. Elegimos un estimador: la media muestral.
  2. Cuantificamos la incertidumbre: el error típico de la media (SE).
  3. Usamos el SE y construimos un intervalo de confianza (IC) para inferir el valor de \(\mu_{\text{edad}}\).
• Un frecuentista concluirá:
  • La media de edad estimada es \(\hat\mu_{\text{edad}}\) con un IC de \([a,b]\).
  • Si repetimos este experimento infinitas veces, el IC contendrá el valor de \(\mu_{\text{edad}}\) en el 95% de los casos.

La edad según el Bayesiano

• La media de edad es una variable aleatoria porque incluimos nuestro conocimiento actual sobre su valor con una distribución de probabilidad.
• El objetivo es combinar (1) los datos y (2) lo que creemos.
  1. Elegimos una distribución a priori para \(\mu_{\text{edad}}\): \(\Pr(\mu_{\text{edad}})\)

La edad según el Bayesiano

• La media de edad es una variable aleatoria porque incluimos nuestro conocimiento actual sobre su valor con una distribución de probabilidad.
• El objetivo es combinar (1) los datos y (2) lo que creemos.
  1. Elegimos una distribución a priori para \(\mu_{\text{edad}}\): \(\Pr(\mu_{\text{edad}})\)
  2. El Teorema de Bayes combina la prior y los datos

La edad según el Bayesiano

• La media de edad es una variable aleatoria porque incluimos nuestro conocimiento actual sobre su valor con una distribución de probabilidad.
• El objetivo es combinar (1) los datos y (2) lo que creemos.
  1. Elegimos una distribución a priori para \(\mu_{\text{edad}}\): \(\Pr(\mu_{\text{edad}})\)
  2. El Teorema de Bayes combina la prior y los datos
  3. Esta combinación es la distribución posterior.

• Un bayesiano concluirá:
  • Tras ver los datos, mi antigua creencia (prior) sobre \(\mu_{\text{edad}}\) se ha actualizado (posterior) .
  • El intervalo de credibilidad al 95% nos dice que hay un 95% de probabilidad de que \(\mu_{\text{edad}}\) esté en \([a,b]\)

Probabilidad como información

Frecuentismo

• La probabilidad de que salga cara es del 50%.

• Esta probabilidad es objetiva:

  • Número de veces que sale cara en infinitas e idénticas repeticiones de un experimento (i.e., el lanzamiento de la moneda).

• La probabilidad es una propiedad “real” del experimento.

Probabilidad como información

Bayesiano

• La probabilidad de que salga cara depende de cuánto sé sobre el proceso de lanzamiento.
• Si no sé nada, diré que ambas cosas son equiprobables (50%)
• A medida que acumulo datos sobre el lanzamiento (ángulo, fuerza, inclinación…), ajusto mi probabilidad acorde al nuevo conocimiento.
• Con información perfecta el resultado está determinado.

Método de estimación

Frecuentismo

• Múltiples métodos:
  • Máxima verosimilitud (clásica, full-information, robusta…)
  • Mínimos cuadrados (ponderados, sin ponderar, robustos…)
• En todos ellos busco un único valor para cada parámetro

Bayesiano

• Estimación basada en simulación: Markov Chain Monte Carlo (MCMC)
• Para cada parámetro obtenemos toda una distribución.
• Suele resumirse esa distribución con la media, la mediana…
• Múltiples algoritmos: Gibbs, NUTS, HMC…

Máxima verosimilitud

Máximos locales y absolutos

Estimación Bayesiana: MCMC

Estimación Bayesiana: MCMC

Estimación Bayesiana: MCMC

Estimación Bayesiana: MCMC

Estimación Bayesiana: MCMC

Estimación Bayesiana: MCMC

Estimación Bayesiana: MCMC

Estimación Bayesiana: MCMC

Estimación Bayesiana: MCMC

Estimación Bayesiana: MCMC

Estimación Bayesiana: MCMC

Estimación Bayesiana: MCMC

Estimación Bayesiana: MCMC

Estimación Bayesiana: MCMC

Estimación Bayesiana: MCMC

• Los algoritmos bayesianos como el sampleador de Gibbs o Monte Carlo Hamiltoniano difieren en cómo seleccionan las muestras.
  • Visualizador de algoritmos MCMC
• En el frecuentismo, el optimizador puede no converger. ¿Cuál es el equivalente a esto en Bayes?
  • Comparación entre varias cadenas de Markov ajustadas a los mismos datos (\(\hat R\)).
  • Calidad de la distribución posterior (tamaño muestral efectivo)
  • Algunos algoritmos, como el Monte Carlo Hamiltoniano, tienen métodos de diagnóstico más avanzados para evaluar la convergencia.

Índices de ajuste:

Frecuentismo

• Índices de ajuste absoluto: \(\chi^2\), RMSEA, SRMR
• Índices de ajuste relativo: CFI, TLI
• Índices de información: AIC, BIC

Bayesiano

• Distribuciones de ajuste absoluto: \(\Gamma\), BRMSEA, BSRMR
• Distribuciones de ajuste relativo: BCFI, BTLI
• Índices de información: WAIC, DIC
• Posterior Predictive Model Checks.

Comparación de modelos

Frecuentismo

• Prueba de razón de verosimilitudes (LRT).
• Incremento en RMSEA y CFI
• Para modelos anidados

Bayesiano

• Comparación de capacidad predictiva: K-fold, LOO
• Comparación de verosimilitudes: factor de Bayes
• Para modelos anidados o no anidados

Modelo teórico

Sintaxis del modelo

• {lavaan} y {blavaan} utilizan la misma sintaxis.
• Especificaremos dos modelos: uno con dos factores y otro con tres factores.

# Sintaxis del modelo: dos factores relacionados
model.2f <- ' visual_text =~ x1 + x2 + x3 + x4 + x5 + x6 
              speed =~ x7 + x8 + x9'
# Sintaxis del modelo: tres factores relacionados
HS.model <- ' visual  =~ x1 + x2 + x3
              textual =~ x4 + x5 + x6
              speed   =~ x7 + x8 + x9 '

# Modelo nulo
null.model <- paste0("x", 1:9, "~~", "x", 1:9, collapse = "\n")

Prior Predictive Check

• ¿Qué priors está utilizando {blavaan} para los parámetros del modelo?

# By-default prior distributions
dpriors()

• Se pueden modificar sus valores, pero no las distribuciones
• Se puede evaluar visualmente la distribución a priori utilizando {blavaan} o R base.

Prior Predictive Check con R base

curve(dnorm(x, mean = 0, sd = 10), from = -50, to = 50, 
      main = "Distribución a priori de Lambda",
      xlab = "Valores de Lambda", ylab = "Densidad de probabilidad", lwd = 2)

Prior Predictive Check con {blavaan}

# Base de datos con variables estandarizadas
df <- HolzingerSwineford1939
df[,paste0("x", 1:9)] <- apply(df[,paste0("x", 1:9)], 2, scale)
# Ajusta un modelo "vacío" prisamp = TRUE!
prior_predictive <- bcfa(model.2f, data = df, std.lv = TRUE, 
                         meanstructure = TRUE, test = "none", sample = 20000, 
                         prisamp = TRUE, bcontrol = list(cores = 3)) 

• Acto seguido, podemos guardar la distribución a priori con la siguiente línea de código

# Guardamos la distribución a priori
prior.draws <- blavInspect(prior_predictive, "mcmc")

• Graficaremos las distribuciones con el paquete {bayesplot}

Prior predictive check: intersecciones

# Prior predictive distribution: intersecciones
mcmc_hist(x = prior.draws, regex_pars = paste0("x", 1:6, "~1"), alpha = .5)

Prior predictive check: pesos factoriales

# Prior predictive distribution: intersecciones
mcmc_hist(x = prior.draws, regex_pars = "visual_text=~", alpha = .5)

Prior predictive check: correlaciones latentes

# Prior predictive distribution: intersecciones
mcmc_hist(x = prior.draws, pars = "visual_text~~speed", alpha = .5)

Prior predictive check: varianzas únicas

# Prior predictive distribution: intersecciones
mcmc_hist(x = prior.draws, pars = paste0("x", 1:6, "~~", "x", 1:6), alpha = .5)

Cambiando la distribución prior

• En {blavaan} no se puede cambiar la distribución a priori, sólo sus parámetros.

# Nueva distribución a priori
new_priors <- dpriors(nu = "normal(0,1)", 
                      lambda = "normal(0,1)", 
                      rho = "beta(2,2)", 
                      psi = "gamma(1,1)")
# Muestreamos de la nueva distribución a priori
new_prior_pred <- bcfa(model.2f, data = df, std.lv = TRUE, meanstructure = TRUE, 
                       test = "none", sample = 20000, prisamp = TRUE,
                       dp = new_priors,
                       bcontrol = list(cores = 3)) 

• Volvemos a extraer la distribución a priori

new.prior.draws <- blavInspect(new_prior_pred, "mcmc")

Prior predictive check: intersecciones

# Prior predictive distribution: intersecciones
mcmc_hist(x = new.prior.draws, regex_pars = paste0("x", 1:6, "~1"), alpha = .5)

Prior predictive check: pesos factoriales

# Prior predictive distribution: intersecciones
mcmc_hist(x = new.prior.draws, regex_pars = "visual_text=~", alpha = .5)

Prior predictive check: correlaciones latentes

# Prior predictive distribution: intersecciones
mcmc_hist(x = new.prior.draws, pars = "visual_text~~speed", alpha = .5)

Prior predictive check: varianzas únicas

# Prior predictive distribution: intersecciones
mcmc_hist(x = new.prior.draws, pars = paste0("x", 1:6, "~~", "x", 1:6), alpha = .5)

Estimación del modelo

• Estimaremos los modelos de dos y tres factores con las nuevas priors.

• Utilizaremos 3 cadenas de Markov con 1500 iteraciones, descartando las primeras 500 iteraciones.

• bcontrol permite paralelizar las cadenas de Markov.

# Estimate all the Bayesian models
blavaan.null.fit <- bcfa(null.model, data = df, burnin = 500, sample = 1000, 
                         meanstructure = TRUE, std.lv = TRUE, 
                         bcontrol = list(cores = 3))
blavaan.2f.fit <- bcfa(model.2f, data = df, burnin = 500, sample = 1000, 
                       meanstructure = TRUE, std.lv = TRUE, 
                       bcontrol = list(cores = 3))
blavaan.HS.fit <- bcfa(HS.model, data = df, burnin = 500, sample = 1000,  
                       meanstructure = TRUE, std.lv = TRUE, 
                       bcontrol = list(cores = 3))

Convergencia del modelo

• La validez de las inferencias depende de la convergencia del modelo.
• Potential Scale Reduction Factor (PSRF, \(\hat R\)) con valores inferiores a 1.05

# Potential Scale Reduction Factor (PSRF): 2 factores
which(blavInspect(blavaan.2f.fit, what = "rhat") > 1.05)

named integer(0)

# Potential Scale Reduction Factor (PSRF): 3 factores
which(blavInspect(blavaan.HS.fit, what = "rhat") > 1.05)

named integer(0)

• Esto sugiere que las cadenas de Markov están muestreando de la misma distribución posterior de manera estable.

Convergencia del modelo

• Podemos evaluar visualmente la convergencia de las cadenas de Markov a través de los gráficos de trayectoria (traceplots).
• Para ello, necesitamos extraer la distribución posterior de los parámetros.

# Distribución posterior de ambos modelos
posterior.2f <- blavInspect(blavaan.2f.fit, "mcmc")
posterior.HS <- blavInspect(blavaan.HS.fit, "mcmc")

• El siguiente código genera traceplots para los pesos factoriales

# Gráficos de trayectoria
mcmc_trace(x = posterior.2f, regex_pars = "visual_text=~")
mcmc_trace(x = posterior.HS, regex_pars = c("visual=~", "textual=~"))

Convergencia del modelo

Convergencia del modelo

Convergencia del modelo

• ¿Qué aspecto tiene este gráfico cuando la convergencia es mala?

Convergencia del modelo

• Podemos destacar una de las cadenas para evaluar su desempeño.

mcmc_trace_highlight(x = posterior.HS, regex_pars = "visual~~textual", highlight = 3)

Eficiencia

• Los valores de los parámetros de iteración a iteración suelen estar correlacionados.
• Cuanto menor sea su correlación, más informativo será el valor en cada muestra.
• Se evalúa a través del tamaño muestral efectivo.
  • Valor mínimo recomendable: \(100\times\) Nº de cadenas de Markov.

min(blavInspect(blavaan.2f.fit, "neff"))

[1] 1427.318

min(blavInspect(blavaan.HS.fit, "neff"))

[1] 1409.651

Ajuste del modelo

• Similares a los índices de ajuste frecuentistas, ¡pero con intervalos de credibilidad!
• Índices Bayesianos de ajuste global
  • Bayesian Root Mean Square Error of Approximation (BRMSEA)
  • Bayesian Unbiased Goodnes-of-fit Index (\(\Gamma\))
  • Bayesian Unbiased Goodnes-of-fit Index ajustado a muestras pequeñas (\(\Gamma_{adj}\))
  • Bayesian McDonald Index (BMc)
• Índices Bayesianos de ajuste incremental
  • Bayesian Comparative Fit Index (BCFI)
  • Bayesian Tucker-Lewis Index (BTLI)

Ajuste del modelo: dos factores

# Indices de ajuste: dos factores
blav_fit_2f <- blavFitIndices(object = blavaan.2f.fit, baseline.model = blavaan.null.fit)
summary(blav_fit_2f, hpd = TRUE, prob = .95, central.tendency = c("mean","median","mode"))

Ajuste del modelo: dos factores

mcmc_hist(data.frame(blav_fit_2f@indices), alpha = 0.5,
          pars = c("BRMSEA", "BCFI", "BTLI", "BGammaHat"))

Ajuste del modelo: tres factores

# Indices de ajuste: tres factores
blav_fit_HS <- blavFitIndices(object = blavaan.HS.fit, baseline.model = blavaan.null.fit)
summary(blav_fit_HS, hpd = TRUE, prob = .95, central.tendency = c("mean","median","mode"))

Ajuste del modelo: tres factores

mcmc_hist(data.frame(blav_fit_HS@indices), alpha = 0.5,
          pars = c("BRMSEA", "BCFI", "BTLI", "BGammaHat"))

Ajuste del modelo: ¿qué hacemos con TLI?

• Como TLI es una distribución, podemos responder preguntas probabilísticas.
• ¿Cuál es la probabilidad de que TLI \(\geq0.90\)?

# Calculamos la función empírica de probabilidad acumulada
TLI_ppost <- ecdf(blav_fit_HS@indices$BTLI)
1 - TLI_ppost(v = 0.90) # Probabilidad de que TLI sea mayor o igual a 0.90

[1] 0.097

# Intervalo de credibilidad al 95%
quantile(blav_fit_HS@indices$BTLI, probs = c(.025, .975))

     2.5%     97.5% 
0.8478141 0.9074505 

• ¿Cuál es la probabilidad de que CFI \(\leq0.90\)?

CFI_ppost <- ecdf(blav_fit_HS@indices$BCFI)
CFI_ppost(v = 0.90) # Probabilidad de que CFI sea menor a 0.90

[1] 0.004

Comparación de modelos

• {blavaan} incluye dos métodos para comparar modelos:

1. Métodos basados en la capacidad predictiva del modelo:
   • Leave-One-Out (LOO) cross-validation.
   • Watanabee-Akaike Information Criteria (WAIC).
   • Ambos incluidos en el paquete {loo}.
   • Cuentan con pruebas de significacion.
2. Métodos basados en la capacidad explicativa del modelo:
   • Factor de Bayes (aproximación de Laplace-Metrópolis).
   • {blavaan} no recomienda su uso.

Comparación de modelos

• La función blavCompare calcula el WAIC y el \(\text{ELPD}_\text{loo}\), una versión más robusta del WAIC.

# Comparación de modelos
blav_com_2f_HS <- blavCompare(object1 = blavaan.2f.fit, object2 = blavaan.HS.fit)
# PSIS-Leave-One-Out: 
blavaan_com_2f_HS$loo[[1]] # modelo de dos factores
blav_com_2f_HS$loo[[2]] # modelo de tres factores

Comparación de modelos

• La diferencia entre los \(\text{ELPD}_\text{loo}\) de ambos modelos, \(\Delta\text{ELPD}_\text{loo}\), permite determinar cuál de los dos tiene mejor capacidad para predecir los datos observados.

# Comparación de ELPDs
blav_com_2f_HS$diff_loo

       elpd_diff se_diff
model2   0.0       0.0  
model1 -45.3      12.0  

• Model2 es el modelo que hemos indicado previamente en object 2

# Comparación de modelos
blav_com_2f_HS <- blavCompare(object1 = blavaan.2f.fit, object2 = blavaan.HS.fit)

• Por tanto, el modelo de tres factores tiene mejor capacidad predictiva que el de dos factores porque la diferencia \(\Delta\mbox{ELPD} = -45.3\) está fuera del rango \(\pm 2 \times 12\).

Parámetros estimados: pesos factoriales

Parámetros estimados: pesos factoriales

mcmc_intervals(posterior.HS, regex_pars = c("visual=~", "textual=~", "speed=~"),  
               prob = 0.8, prob_outer = 0.95, point_est = "median")

Parámetros estimados: correlaciones latentes

Parámetros estimados: correlaciones latentes

mcmc_areas(posterior.HS, regex_pars = c("visual~~", "textual~~", "speed~~"),  
           prob = 0.8, prob_outer = 0.95, point_est = "median")

Posterior Predictive Model Checks (PPMC)

• Método para evaluar la capacidad predictiva de un modelo.
• La idea central es verificar si el model puede generar datos similares a los datos observados.
• ¿Podemos distinguir los datos observados de datos simulados desde el modelo?
• En {blavaan} podemos utilizar la función sampleData para generar datos desde el modelo.
• Tras simularlo, se comparan gráficamente los datos observados y simulados
• Finalmente, podemos ver el parecido entre \(y\) con \(y^{\text{rep}}\)

Posterior Predictive Model Checks: \(y^{\text{rep}}\)

# Simulamos 100 bases de datos
yrep <- sampleData(bcfa.HS.fit, nrep = 100, simplify = TRUE)
# Guardamos las respuestas de los ítems en un data frame a parte
Y <- df[,7:15]
# Guardamos los resultados para cada item
posterior_plots_list <- vector(mode = "list", length = ncol(Y))
for(i in 1:ncol(Y)){
  posterior_plots_list[[i]] <- ppc_dens_overlay(
    y = Y[,i], yrep = t(sapply(yrep, function(x) x[,i]))
    ) + 
    labs(title = paste("Item", i))
}
# Gráfico comparando y e yrep
cowplot::plot_grid(plotlist = posterior_plots_list)

Posterior Predictive Model Checks: \(y^{\text{rep}}\)

Posterior Predictive Model Checks (PPMC)

Una iteración

• Tenemos \(\color{#0197FD}{\boldsymbol{\nu}^{(1)}}\) y \(\color{#0197FD}{\mathbf{\Sigma^{(1)}}} = \mathbf{\Lambda^{(1)}}\cdot\mathbf{\Psi^{(1)}}\cdot\mathbf{\Lambda^{\prime(1)}} + \mathbf{\Theta^{(1)}}\)
• Calculamos \(\chi^{2(1)}\): \(\text{Cov}(\mathbf{Y}) - \color{#0197FD}{\mathbf{\Sigma^{(1)}}}\)
• Simulamos datos utilizando estos parámetros: \(\mathbf{Y}_{\text{sim}}^{(1)}\)

• Calculamos \(\chi^{2(1)}_\text{sim}\): \(\text{Cov}\left(\mathbf{Y}_{\text{sim}}^{(1)}\right) - \color{#0197FD}{\mathbf{\Sigma^{(1)}}}\)

• Comparamos la distribución de \(\chi^{2(1)}\) y \(\chi^{2(1)}_\text{sim}\).

Posterior Predictive Model Checks (PPMC)

Posterior Predictive Model Checks (PPMC)

Posterior Predictive Model Checks (PPMC)

• El PPMC puede aplicarse a cualquier estadístico.
• Durante la estimación, en cada iteración:
  1. Se extrae una muestra de los parámetros \(\theta^{(s)}\) de la distribución posterior (esto lo hacemos siempre).
  2. Se utilizan estos parámetros para generar datos simulados: \(\mathbf{Y}_{\text{sim}}\).
  3. Se calcula algún estadístico de interés en \(\mathbf{Y}\) e \(\mathbf{Y}_{\text{sim}}\).
  4. Se compara la distribución del estadístico entre \(\mathbf{Y}\) e \(\mathbf{Y}_{\text{sim}}\).
• En {blavaan} pueden realizarse con la función ppmc

# Posterior Predictive Model Checks
ppmc_fit_2f <- ppmc(object = bcfa.2f.fit, thin = 1, fit.measures = c("srmr","chisq"))
ppmc_fit_HS <- ppmc(object = bcfa.HS.fit, thin = 1, fit.measures = c("srmr","chisq"))

Posterior Predictive Model Checks (PPMC): χ²

Posterior Predictive Model Checks (PPMC): χ²

Posterior Predictive Model Checks (PPMC): χ²

Posterior Predictive Model Checks (PPMC): χ²

Posterior Predictive Model Checks (PPMC): SRMR

Posterior Predictive Model Checks (PPMC): SRMR

Posterior Predictive Model Checks (PPMC): SRMR

Posterior Predictive Model Checks (PPMC): SRMR

PPMC y Fiabilidad: ω

• Podemos utilizar ppmc con cualquier función que pueda utilizarse para un modelo de {lavaan}.
• Podemos utilizar compRelSEM de {SemTools} para obtener la distribución posterior de \(\omega\).

# Función para obtener la fiabilidad en cada iteración
bcfa.reliability <- function(fit){ semTools::compRelSEM(fit) }
# Fiabilidad en el modelo unidimensional
ppmc_omega_2f <- ppmc(object = blavaan.2f.fit, discFUN = bcfa.reliability)
# Fiabilidad en el modelo de dos factores
ppmc_omega_HS <- ppmc(object = blavaan.HS.fit, discFUN = bcfa.reliability)
# La distribución posterior se guarda aquí
posterior_omega_HS@obsDist$discFUN1

PPMC y Fiabilidad: ω

PPMC y Fiabilidad: ω

Psicometría Bayesiana

Únete al lado oscuro de la fuerza

1. Está emergiendo con más fuerza que nunca.
2. Es un marco de trabajo coherente.
3. Ideal para muestras pequeñas (i.e., \(N=70\)).
4. Incertidumbre en cualquier parámetro.
5. En última instancia, es frecuentismo liberado.
6. Puedes especificar y ajustar cualquier modelo que puedas imaginar.